{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 神经元网络\n",
    "\n",
    "## 1.BP神经网络\n",
    "见 'MPLClassifier.ipynb'\n",
    "\n",
    "## 2.RBF神经网络\n",
    "\n",
    "- RBF神经网络经过这样两层变化: $\\left\\{\\begin{array}{l}R_i(X)=exp(-||X-C_i||^2/2\\sigma_i^2),\\qquad i=1,\\cdots,m\\\\\\hat{y}_k=\\sum\\limits_{i=1}^{m}\\omega_{ik}R_i(X),\\qquad k=1,\\cdots,p\\end{array}\\right.$  \n",
    "这样只有小部分靠近中心的隐藏层神经元被激活($R_i(X)$随着其中范数增大,指数减少)\n",
    "- 确定基函数中心$C_i$. 一般采用K均值聚类法.\n",
    "- 确定基函数宽度$\\sigma_i$. 通常令它等于基函数中心与子样本集中样本模式之间的平均距离\n",
    "- 确定权值$\\omega_{ik}$. 采用最小均方误差测度."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "WARNING:tensorflow:From D:\\Anaconda3\\lib\\site-packages\\tensorflow_core\\python\\ops\\resource_variable_ops.py:1635: calling BaseResourceVariable.__init__ (from tensorflow.python.ops.resource_variable_ops) with constraint is deprecated and will be removed in a future version.\n",
      "Instructions for updating:\n",
      "If using Keras pass *_constraint arguments to layers.\n",
      "Loss function at step 0 is 456.01239013671875\n",
      "Loss function at step 20 is 4.353085041046143\n",
      "Loss function at step 40 is 2.0776684284210205\n",
      "Loss function at step 60 is 1.0662990808486938\n",
      "Loss function at step 80 is 0.7423067688941956\n",
      "Loss function at step 100 is 0.6860532164573669\n",
      "Loss function at step 120 is 0.7009456157684326\n",
      "Loss function at step 140 is 0.7049466967582703\n",
      "Loss function at step 160 is 0.6740176677703857\n",
      "Loss function at step 180 is 0.6010324358940125\n",
      "Training complete!\n",
      "Relative Error Sum: 134.1192%\n",
      "                    Input Pred Out Raw Out     Error Relative Error\n",
      "0   [15.6  5.6  3.5 25.5]  22.8648    22.9   -0.0352       -0.1536%\n",
      "1   [27.8  4.3  1.   7.7]  23.3983    23.4   -0.0017       -0.0071%\n",
      "2   [35.2  3.  38.1  3.7]  36.6834    36.8   -0.1166       -0.3168%\n",
      "3   [10.2  3.4  3.5  7.4]  20.2288    22.0   -1.7712       -8.0509%\n",
      "4   [29.1 33.2  1.6 24. ]   6.4763     6.4    0.0763        1.1918%\n",
      "5   [10.2 11.6  2.2 26.7]  29.4712    29.4    0.0712        0.2423%\n",
      "6   [35.4  4.1  1.3  7. ]  24.4484    26.2   -1.7516       -6.6856%\n",
      "7       [8.7 3.5 7.5 5. ]  20.7974    20.9   -0.1026       -0.4909%\n",
      "8   [25.4  0.7 22.2 35.4]  26.4857    26.5   -0.0143       -0.0541%\n",
      "9   [15.3  6.   2.  17.5]  37.3170    37.3    0.0170        0.0456%\n",
      "10  [25.9  1.2  9.   3.3]  25.1844    22.8    2.3844       10.4580%\n",
      "11  [64.3  3.7  4.6  4.8]  18.6040    19.8   -1.1960       -6.0402%\n",
      "12  [55.9  2.9  0.3  5.2]  20.8820    19.6    1.2820        6.5408%\n",
      "13  [19.6 10.5 10.7 10.3]  23.9080    28.5   -4.5920      -16.1124%\n",
      "14  [35.6  2.4  6.6 24.6]  22.7483    22.8   -0.0517       -0.2268%\n",
      "15  [10.9  9.4  0.8  7.1]  19.9576    18.2    1.7576        9.6569%\n",
      "16  [24.7  8.2  7.7 14.4]  23.2665    23.8   -0.5335       -2.2416%\n",
      "17  [22.6 11.2  9.9 18.5]  21.7538    17.3    4.4538       25.7445%\n",
      "18  [21.5  2.9  1.6  4.5]  22.6010    21.9    0.7010        3.2010%\n",
      "19  [54.7  3.3  3.7 11.6]  20.7761    32.8  -12.0239      -36.6581%\n"
     ]
    }
   ],
   "source": [
    "# RBF神经网络(based on tensorflow): https://github.com/shiluqiang/RBF_NN_tensorflow/blob/master/RBF_tensorflow.py\n",
    "# 由于数据量小, 本人的不熟悉等原因, 预测结果并不理想.\n",
    "\n",
    "import tensorflow as tf\n",
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "import numpy as np\n",
    "from mpl_toolkits.mplot3d import Axes3D\n",
    "from sklearn.cluster import KMeans\n",
    "import pandas as pd\n",
    "from sklearn import preprocessing\n",
    "\n",
    "class RBF_NN():\n",
    "    def __init__(self, hidden_nodes, input_data_trainX, input_data_trainY):\n",
    "        self.hidden_nodes = hidden_nodes #隐含层节点数\n",
    "        self.input_data_trainX = input_data_trainX #训练样本的特征\n",
    "        self.input_data_trainY = input_data_trainY #训练样本的标签\n",
    "    \n",
    "    def fit(self):\n",
    "        '''模型训练\n",
    "        '''\n",
    "        tf.compat.v1.disable_eager_execution()\n",
    "        # 1.声明输入输出的占位符\n",
    "        n = 19\n",
    "        n_input = (self.input_data_trainX).shape[1]\n",
    "        n_output = (self.input_data_trainY).shape[0]\n",
    "        X = tf.compat.v1.placeholder('float', [None, n_input],name = 'X')\n",
    "        Y = tf.compat.v1.placeholder('float', [None, 1],name = 'Y')\n",
    "        \n",
    "        # 2.参数设置\n",
    "        ## RBF函数参数\n",
    "        ### K-means求中心\n",
    "        random_state = 170\n",
    "        kms = KMeans(n_clusters=self.hidden_nodes, random_state=None)\n",
    "        pred = kms.fit_predict(trainX)\n",
    "        \n",
    "        # c = tf.Variable(tf.random_normal(shape=(self.hidden_nodes, n_input)),name = 'c')\n",
    "        # c = tf.concat((tf.cast(tf.Variable(kms.cluster_centers_), tf.float32), c), axis=0)\n",
    "        c = tf.cast(tf.Variable(kms.cluster_centers_), tf.float32)\n",
    "        delta = tf.Variable(tf.compat.v1.random_normal(shape=(1,self.hidden_nodes)), name='delta')\n",
    "        ## 隐含层到输出层权重和偏置\n",
    "        W = tf.Variable(tf.compat.v1.random_normal(shape=(self.hidden_nodes, 1)), name='W')\n",
    "        b = tf.Variable(tf.compat.v1.random_normal(shape=(1, 1)), name='b')\n",
    "        \n",
    "        # 3.构造前向传播计算图\n",
    "        ## 隐含层输出\n",
    "        ### 特征样本与RBF均值的距离\n",
    "        dist = tf.reduce_sum(tf.square(tf.subtract(tf.tile(X,[self.hidden_nodes, 1]),c)), axis=1)\n",
    "        dist = tf.multiply(1.0,tf.transpose(dist))\n",
    "        ### RBF方差的平方\n",
    "        delta_2 = tf.square(delta)\n",
    "        ### 隐含层输出\n",
    "        RBF_OUT = tf.exp(tf.multiply(-1.0,tf.divide(dist,tf.multiply(2.0,delta_2))))\n",
    "        ## 输出层输入\n",
    "        output_in = tf.matmul(RBF_OUT, W) + b\n",
    "        \n",
    "        # 4.声明代价函数优化算法\n",
    "        loss = tf.reduce_mean(tf.pow(Y - output_in,2)) #损失函数为均方误差\n",
    "        train_op = tf.compat.v1.train.AdamOptimizer(0.1).minimize(loss) #优化算法为梯度下降法\n",
    "    \n",
    "        # 5.反向传播求参数\n",
    "        trX = self.input_data_trainX[:n]\n",
    "        trY = self.input_data_trainY[:n]\n",
    "        \n",
    "        with tf.compat.v1.Session() as sess:\n",
    "            ## 初始化所有参数\n",
    "            tf.compat.v1.global_variables_initializer().run()\n",
    "            for epoch in range(200):\n",
    "                for i in range(trX.shape[0]):\n",
    "                    feed = {X:trX[i][:,None].T, Y:[[trY[i]]]}\n",
    "                    sess.run(train_op,feed_dict=feed)\n",
    "                if epoch % 20. == 0 :\n",
    "                    total_loss = 0.0\n",
    "                    for j in range(trX.shape[0]):\n",
    "                        total_loss += sess.run(loss, feed_dict={X:trX[i][:,None].T, Y:[[trY[i]]]})\n",
    "                    print('Loss function at step %d is %s'%(epoch, total_loss / trX.shape[0]))\n",
    "                    \n",
    "            print('Training complete!')\n",
    "\n",
    "            W = W.eval()\n",
    "            b = b.eval()\n",
    "            c = c.eval()\n",
    "            delta = delta.eval()\n",
    "            pred_trX = np.mat(np.zeros((len(trX),n_output)))\n",
    "            \n",
    "            ## 训练准确率\n",
    "            correct_tr = 0.0\n",
    "            pred = []\n",
    "            for i in range(self.input_data_trainX.shape[0]):\n",
    "                pred_tr = sess.run(output_in, feed_dict={X:self.input_data_trainX[i][:,None].T})\n",
    "                pred.append(pred_tr[0][0])\n",
    "            df_columns = ['Input', 'Pred Out', 'Raw Out', 'Error', 'Relative Error']\n",
    "            pred = np.array(pred)\n",
    "            df = pd.DataFrame(np.c_[[x.__str__() for x in self.input_data_trainX],\n",
    "                                    pred,\n",
    "                                    self.input_data_trainY,\n",
    "                                    np.subtract(pred, self.input_data_trainY),\n",
    "                                    np.array(np.divide(np.subtract(pred, self.input_data_trainY), self.input_data_trainY), dtype=np.float)\n",
    "                                   ], columns=df_columns)\n",
    "            print('Relative Error Sum: {:.4%}'.format(df['Relative Error'].astype(float).abs().sum()))\n",
    "            df['Pred Out'] = df['Pred Out'].apply(lambda x: format(float(x), '.4f'))\n",
    "            df['Raw Out'] = df['Raw Out'].apply(lambda x: format(float(x), '.1f'))\n",
    "            df['Error'] = df['Error'].apply(lambda x: format(float(x), '.4f'))\n",
    "            df['Relative Error'] = df['Relative Error'].apply(lambda x: format(float(x), '.4%'))\n",
    "            print(df)\n",
    "\n",
    "data = np.loadtxt('15.D 水库年径流与因子特征.txt')\n",
    "# data = preprocessing.scale(data)\n",
    "\n",
    "trainX = data[:, :-1]\n",
    "trainY = data[:, -1]\n",
    "input_data_trainX = trainX\n",
    "input_data_trainY = trainY\n",
    "rbf = RBF_NN(10, input_data_trainX, input_data_trainY)\n",
    "rbf.fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "At the end of epoch: 0, with loss: 1282.6154378255208\n",
      "At the end of epoch: 200, with loss: 46.1252498626709\n",
      "At the end of epoch: 400, with loss: 28.372927983601887\n",
      "At the end of epoch: 600, with loss: 17.97803243001302\n",
      "At the end of epoch: 800, with loss: 10.040056864420572\n",
      "Relative Error Sum: 973.1400%\n",
      "                    Input Pred Out Raw Out    Error Relative Error\n",
      "0   [15.6  5.6  3.5 25.5]  36.6835    22.9  13.7835       60.1899%\n",
      "1   [27.8  4.3  1.   7.7]  22.5843    23.4  -0.8157       -3.4860%\n",
      "2   [35.2  3.  38.1  3.7]  54.3693    36.8  17.5693       47.7428%\n",
      "3   [10.2  3.4  3.5  7.4]  12.0602    22.0  -9.9398      -45.1808%\n",
      "4   [29.1 33.2  1.6 24. ]  52.4762     6.4  46.0762      719.9399%\n",
      "5   [10.2 11.6  2.2 26.7]  31.6173    29.4   2.2173        7.5419%\n",
      "6   [35.4  4.1  1.3  7. ]  26.2944    26.2   0.0944        0.3604%\n",
      "7       [8.7 3.5 7.5 5. ]  19.0089    20.9  -1.8911       -9.0485%\n",
      "8   [25.4  0.7 22.2 35.4]  26.5318    26.5   0.0318        0.1199%\n",
      "9   [15.3  6.   2.  17.5]  31.2433    37.3  -6.0567      -16.2378%\n",
      "10  [25.9  1.2  9.   3.3]  22.2908    22.8  -0.5092       -2.2335%\n",
      "11  [64.3  3.7  4.6  4.8]  21.8162    19.8   2.0162       10.1829%\n",
      "12  [55.9  2.9  0.3  5.2]  17.3938    19.6  -2.2062      -11.2560%\n",
      "13  [19.6 10.5 10.7 10.3]  27.3850    28.5  -1.1150       -3.9124%\n",
      "14  [35.6  2.4  6.6 24.6]  23.2914    22.8   0.4914        2.1553%\n",
      "15  [10.9  9.4  0.8  7.1]  19.3219    18.2   1.1219        6.1643%\n",
      "16  [24.7  8.2  7.7 14.4]  23.7118    23.8  -0.0882       -0.3707%\n",
      "17  [22.6 11.2  9.9 18.5]  20.7248    17.3   3.4248       19.7968%\n",
      "18  [21.5  2.9  1.6  4.5]  20.8946    21.9  -1.0054       -4.5909%\n",
      "19  [54.7  3.3  3.7 11.6]  33.6625    32.8   0.8625        2.6295%\n"
     ]
    }
   ],
   "source": [
    "# BP神经网络(based on keras)\n",
    "\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "from tensorflow.keras.models import Sequential\n",
    "from tensorflow.keras.layers import Dense, Activation\n",
    "from tensorflow.keras.callbacks import LambdaCallback\n",
    "\n",
    "data = np.loadtxt('15.D 水库年径流与因子特征.txt')\n",
    "X = data[:, :-1]\n",
    "Y = data[:, -1]\n",
    "n = 5\n",
    "trainX = X[n:]\n",
    "trainY = Y[n:]\n",
    "\n",
    "model = Sequential()  #层次模型\n",
    "model.add(Dense(12, input_dim=4)) #输入层，Dense表示BP层\n",
    "model.add(Activation('relu'))\n",
    "model.add(Dense(5, input_dim=12))\n",
    "model.add(Activation('relu'))\n",
    "model.add(Dense(1, input_dim=5))  #输出层\n",
    "model.compile(loss='mean_squared_error', optimizer='Adam') #编译模型\n",
    "\n",
    "\n",
    "def eposh_callback(epoch, logs):\n",
    "    if epoch % 200 == 0:\n",
    "        print('At the end of epoch: {}, with loss: {}'.format(epoch, logs['loss']))\n",
    "        \n",
    "        \n",
    "batch_print_callback = LambdaCallback(on_epoch_end=eposh_callback)\n",
    "model.fit(trainX, trainY, epochs = 1000, batch_size = 5, verbose=0, callbacks=[batch_print_callback]) #训练模型1000次\n",
    "\n",
    "pred = np.array(model.predict(X)).flatten()\n",
    "df_columns = ['Input', 'Pred Out', 'Raw Out', 'Error', 'Relative Error']\n",
    "df = pd.DataFrame(np.c_[[x.__str__() for x in X],\n",
    "                    pred,\n",
    "                    Y,\n",
    "                    np.subtract(pred, Y),\n",
    "                    np.array(np.divide(np.subtract(pred, Y), Y), dtype=np.float)\n",
    "                   ], columns=df_columns)\n",
    "print('Relative Error Sum: {:.4%}'.format(df['Relative Error'].astype(float).abs().sum()))\n",
    "df['Pred Out'] = df['Pred Out'].apply(lambda x: format(float(x), '.4f'))\n",
    "df['Raw Out'] = df['Raw Out'].apply(lambda x: format(float(x), '.1f'))\n",
    "df['Error'] = df['Error'].apply(lambda x: format(float(x), '.4f'))\n",
    "df['Relative Error'] = df['Relative Error'].apply(lambda x: format(float(x), '.4%'))\n",
    "print(df)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
